The data contained in these files are described in detail in,
https://www.science.org/doi/10.1126/scitranslmed.abd8995
For convenience we are providing 3 processed data tables here. One is a raw sparse matrix, that is the output of our processing pipeline with no further manipulation. We are also providing a fully processed UMI table that has been filtered, dimension reduced, clustered, and normalized. The processed UMI table is provided in the expression set class format (https://www.bioconductor.org/packages/devel/bioc/vignettes/Biobase/inst/doc/ExpressionSetIntroduction.pdf). We are also providing a separate processed UMI table for the T Cells.
The raw .fastq files have been deposited through dbGAP.
https://www.ncbi.nlm.nih.gov/projects/gap/cgi-bin/study.cgi?study_id=phs002455.v1.p1
These files are protected, and you will need an account and to request access to download the files. As of 10/5/21 the data files are processing on the dbGAP server.
# devtools::install_github("garber-lab/SignallingSingleCell")
library(SignallingSingleCell)
# load(url("https://www.dropbox.com/s/qkk4kzsaoz3seoj/phs002455_raw_UMITable.Rdata?dl=1")) # The Raw Data
load(url("https://www.dropbox.com/s/ipr7dys1ixd7f7m/phs002455_processed_UMItable.Rdata?dl=1"))
load(url("https://www.dropbox.com/s/yofmvzfb36j64lv/phs002455_TCells_processed_UMItable.Rdata?dl=1"))
The ExpressionSet class (ex_sc) is a convenient data structure that contains 3 dataframes. These dataframes contain expression data, cell information, and gene information respectively.
exprs(ex_sc) is the expression data, where rows are genes and columns are cells.
pData(ex_sc) is cell information, where rows are cells and columns are metadata.
fData(ex_sc) is gene information, where rows are genes and columns are metadata.
colnames(pData(phs002455_processed_UMItable))
[1] "Patient" "Disease" "Phenotype" "Skin"
[5] "UMI_sum_raw" "x" "y" "Cluster"
[9] "UMI_sum_norm" "CellType" "Cluster_Fig_1B" "Cluster_Refined"
# Columns 1-4 are all phenotypic metadata related to the patient sample
# Columns 5-12 are all related to cell specific metadata such as tSNE coordinates, Clusters, etc
colnames(fData(phs002455_processed_UMItable))
[1] "CD8_marker_score_Cluster_Refined"
[2] "DC_marker_score_Cluster_Refined"
[3] "GD_marker_score_Cluster_Refined"
[4] "KRT-B1_marker_score_Cluster_Refined"
[5] "KRT-B2_marker_score_Cluster_Refined"
[6] "KRT-ECR_marker_score_Cluster_Refined"
[7] "KRT-GR_marker_score_Cluster_Refined"
[8] "KRT-SP_marker_score_Cluster_Refined"
[9] "MAC_marker_score_Cluster_Refined"
[10] "MEL_marker_score_Cluster_Refined"
[11] "NK_marker_score_Cluster_Refined"
[12] "TConv_marker_score_Cluster_Refined"
[13] "Treg_marker_score_Cluster_Refined"
# The only thing stored in fData right now are cluster marker scores.
# These were calculated with id_markers()
Reproduce Figure 1B from https://www.science.org/doi/10.1126/scitranslmed.abd8995
plot_tsne_metadata(phs002455_processed_UMItable, color_by = "Cluster_Fig_1B", shuffle = T)
plot_tsne_gene(phs002455_processed_UMItable, gene = c("TRAC", "TYR"))
# Additional plots
plot_violin(phs002455_processed_UMItable, gene = c("IFNG"), color_by = "Skin", facet_by = "Cluster_Refined")
markers <- return_markers(phs002455_processed_UMItable, return_by = "Cluster_Refined", num_markers = 5)
markers
$CD8_Markers
[1] "CD8A" "NKG7" "CCL5" "CD8B" "CD3D"
$DC_Markers
[1] "CD207" "CD1A" "FCER1A" "HLA-DQB2" "AIF1"
$GD_Markers
[1] "TARP" "TRGC1" "CD8A" "CCL5" "TRGC2"
$`KRT-B1_Markers`
[1] "CYR61" "COL17A1" "POSTN" "SYT8" "RND3"
$`KRT-B2_Markers`
[1] "S100A2" "KRT15" "CCDC3" "TP63" "MIR205HG"
$`KRT-ECR_Markers`
[1] "KRT6B" "MUCL1" "S100P" "MMP7" "KRT77"
$`KRT-GR_Markers`
[1] "KRT2" "LGALS7" "SBSN" "KLK11" "KRTDAP"
$`KRT-SP_Markers`
[1] "KRTDAP" "LGALS7" "LGALS7B" "SBSN" "DMKN"
$MAC_Markers
[1] "CLEC10A" "FPR3" "LGALS2" "LYZ" "CPVL"
$MEL_Markers
[1] "APOD" "TYR" "BCAN" "DCT" "PLP1"
$NK_Markers
[1] "SPINK2" "XCL1" "CTSW" "KLRB1" "TRDC"
$TConv_Markers
[1] "IL32" "TRBC2" "CD3D" "CD6" "CD48"
$Treg_Markers
[1] "CD27" "CD3D" "IL32" "TRBC2" "CTLA4"
plot_gene_dots(phs002455_processed_UMItable, genes = unique(unlist(markers)), break_by = "Cluster_Refined")
Reproduce Figure 1C from https://www.science.org/doi/10.1126/scitranslmed.abd8995
plot_tsne_metadata(phs002455_TCells_processed_UMItable, color_by = "Cluster", shuffle = T)
plot_tsne_gene(phs002455_TCells_processed_UMItable, gene = c("TRAC", "CD4", "CD8A", "FOXP3", "TRGC1", "FCER1G"))
# To create the aggregate bulk heatmaps first calculate the aggregate bulk values
phs002455_TCells_processed_UMItable <- calc_agg_bulk(phs002455_TCells_processed_UMItable, aggregate_by = "Cluster")
# Then identify markers
phs002455_TCells_processed_UMItable <- id_markers(phs002455_TCells_processed_UMItable, id_by = "Cluster", overwrite = T)
[1] "Finding markers based on fraction expressing"
[1] "Finding markers based on mean expressing"
[1] "Merging Lists"
markers <- return_markers(phs002455_TCells_processed_UMItable)
markers
$CD8_Markers
[1] "CD8A" "GZMH" "CD8B" "NKG7" "GZMK" "GZMB" "EOMES" "GZMA" "CCL4"
[10] "CCL5"
$GD_Markers
[1] "TARP" "TRGC1" "TRGC2" "KLRC3" "CD8A" "CD8B" "NKG7"
[8] "FXYD2" "KIR3DL2" "KLRC2"
$NK_Markers
[1] "SPINK2" "FCER1G" "XCL2" "KLRF2" "XCL1" "MB"
[7] "TMPRSS11E" "KRT81" "SPTSSB" "TRDC"
$TConv_Markers
[1] "CD40LG" "CD4" "CD6" "CCR6" "TNFRSF25" "ADAM19"
[7] "IL26" "CTSH" "HLF" "RORA"
$Treg_Markers
[1] "CD4" "FOXP3" "CTLA4" "GBP5" "TBC1D4" "TIGIT" "CD27" "LAIR2"
[9] "IKZF4" "F5"
markers <- unique(unlist(markers))
Now create a heatmap using these values
plot_heatmap(phs002455_TCells_processed_UMItable, type = "bulk", genes = markers)
These functions are still in development, but we have included them for illustration purposes. The goal is to take aggregate bulk CPM values derived from scRNA-seq, cross reference them to a database of ligand and receptor pairs (Ramilowski, Jordan A et al. “A draft network of ligand-receptor-mediated multicellular signalling in human.” Nature communications vol. 6 7866. 22 Jul. 2015, doi:10.1038/ncomms8866), and then construct a network. The first step is to calculate the aggregate bulk CPM values, followed by annotating the ligands and receptors in the data.
# This function is calculating our aggregate bulk CPM values. We have some thresholds (25 CPM or the value will be set to 0, and 15 minimum cells expressing within the group or the value set to 0). This is to greatly reduce the complexity of the network, and reduce spurious nodes with little supporting data.
phs002455_processed_UMItable <- calc_agg_bulk(phs002455_processed_UMItable, aggregate_by = c("Skin", "Cluster_Refined"), group_by = "Skin", cutoff_cpm = 25, cutoff_num = 15)
# The aggregate bulk data is written into fData()
colnames(fData(phs002455_processed_UMItable))[14:ncol(fData(phs002455_processed_UMItable))]
[1] "1-Healthy_CD8_num_genes_11521_num_cells_311_percent_3.28_bulk"
[2] "2-NonLesional_CD8_num_genes_10925_num_cells_197_percent_1.79_bulk"
[3] "3-Disease_CD8_num_genes_13319_num_cells_1444_percent_12.13_bulk"
[4] "1-Healthy_DC_num_genes_12950_num_cells_1860_percent_19.64_bulk"
[5] "2-NonLesional_DC_num_genes_11836_num_cells_753_percent_6.82_bulk"
[6] "3-Disease_DC_num_genes_12558_num_cells_1184_percent_9.95_bulk"
[7] "1-Healthy_GD_num_genes_9552_num_cells_96_percent_1.01_bulk"
[8] "2-NonLesional_GD_num_genes_9891_num_cells_116_percent_1.05_bulk"
[9] "3-Disease_GD_num_genes_12381_num_cells_564_percent_4.74_bulk"
[10] "1-Healthy_KRT-B1_num_genes_8841_num_cells_157_percent_1.66_bulk"
[11] "2-NonLesional_KRT-B1_num_genes_8561_num_cells_142_percent_1.29_bulk"
[12] "3-Disease_KRT-B1_num_genes_10920_num_cells_401_percent_3.37_bulk"
[13] "1-Healthy_KRT-B2_num_genes_12948_num_cells_1909_percent_20.16_bulk"
[14] "2-NonLesional_KRT-B2_num_genes_13707_num_cells_3924_percent_35.57_bulk"
[15] "3-Disease_KRT-B2_num_genes_13762_num_cells_3721_percent_31.27_bulk"
[16] "1-Healthy_KRT-ECR_num_genes_7783_num_cells_108_percent_1.14_bulk"
[17] "2-NonLesional_KRT-ECR_num_genes_6756_num_cells_78_percent_0.71_bulk"
[18] "3-Disease_KRT-ECR_num_genes_8872_num_cells_171_percent_1.44_bulk"
[19] "1-Healthy_KRT-GR_num_genes_12351_num_cells_1100_percent_11.61_bulk"
[20] "2-NonLesional_KRT-GR_num_genes_11157_num_cells_457_percent_4.14_bulk"
[21] "3-Disease_KRT-GR_num_genes_10234_num_cells_193_percent_1.62_bulk"
[22] "1-Healthy_KRT-SP_num_genes_13079_num_cells_2593_percent_27.38_bulk"
[23] "2-NonLesional_KRT-SP_num_genes_13648_num_cells_4431_percent_40.16_bulk"
[24] "3-Disease_KRT-SP_num_genes_12802_num_cells_1871_percent_15.72_bulk"
[25] "1-Healthy_MAC_num_genes_10720_num_cells_75_percent_0.79_bulk"
[26] "2-NonLesional_MAC_num_genes_10462_num_cells_63_percent_0.57_bulk"
[27] "3-Disease_MAC_num_genes_12454_num_cells_192_percent_1.61_bulk"
[28] "1-Healthy_MEL_num_genes_13299_num_cells_591_percent_6.24_bulk"
[29] "2-NonLesional_MEL_num_genes_13235_num_cells_503_percent_4.56_bulk"
[30] "3-Disease_MEL_num_genes_11742_num_cells_151_percent_1.27_bulk"
[31] "1-Healthy_NK_num_genes_11275_num_cells_228_percent_2.41_bulk"
[32] "2-NonLesional_NK_num_genes_9811_num_cells_80_percent_0.73_bulk"
[33] "3-Disease_NK_num_genes_12048_num_cells_365_percent_3.07_bulk"
[34] "1-Healthy_TConv_num_genes_11591_num_cells_303_percent_3.2_bulk"
[35] "2-NonLesional_TConv_num_genes_10372_num_cells_145_percent_1.31_bulk"
[36] "3-Disease_TConv_num_genes_12256_num_cells_571_percent_4.8_bulk"
[37] "1-Healthy_Treg_num_genes_10566_num_cells_140_percent_1.48_bulk"
[38] "2-NonLesional_Treg_num_genes_10429_num_cells_144_percent_1.31_bulk"
[39] "3-Disease_Treg_num_genes_13003_num_cells_1073_percent_9.02_bulk"
# This function will now cross reference the ligand and receptor database, and then annotate this information into fData
phs002455_processed_UMItable <- id_rl(phs002455_processed_UMItable)
# The ligand and receptor annotations are written into fData(). The IFNG entries are highlighted below.
fData(phs002455_processed_UMItable)[c("IFNG", "IFNGR1", "IFNGR2"),c(53:56)]
networks_Receptors networks_ligands networks_ligs_to_receptor
IFNG FALSE TRUE <NA>
IFNGR1 TRUE FALSE IFNG
IFNGR2 TRUE FALSE IFNG
networks_receptor_to_ligs
IFNG IFNGR1_IFNGR2
IFNGR1 <NA>
IFNGR2 <NA>
vals <- which(fData(phs002455_processed_UMItable)[,"networks_ligands"])
ligs <- rownames(fData(phs002455_processed_UMItable))[vals]
length(ligs)
[1] 278
plot_heatmap(phs002455_processed_UMItable, genes = ligs, type = "bulk", cluster_by = "row",facet_by = "Cluster_Refined", pdf_format = "tile", scale_by = "row", cluster_type = "kmeans", k = 15)
vals <- which(fData(phs002455_processed_UMItable)[,"networks_Receptors"])
recs <- rownames(fData(phs002455_processed_UMItable))[vals]
length(recs)
[1] 319
plot_heatmap(phs002455_processed_UMItable, genes = recs, type = "bulk", cluster_by = "row",facet_by = "Cluster_Refined", pdf_format = "tile", scale_by = "row", cluster_type = "kmeans", k = 15)
Now we can calculate the long format network table. Please note this step is slow and poorly optimized… it needs to be rewritten. On my laptop it takes about 30 minutes to run.
Please keep in mind the network scales exponentially with the number of groups you provide. For each individual ligand or receptor, a node is created when the cell type expresses it, and an edge between all of its expressed pairs. This means that if 5 cells express ligand A which binds to receptor B, there will be 5 ligand A nodes, and an edge from each of them to receptor B node. The more promiscuous a given ligand and receptor interaction, and the more cell types in the data, this network quickly explodes to become unmanageable in size.
We have found 10-15 cell types is the upper limit of what can be done on a local computer without needing a cluster environment.
For convenience a download link to the output file is provided below, however you may opt to run the command on your own.
### Warning Slow!!! Use provided output table for speed. The command is left commented for speed reasons ###
# network_table <- calc_rl_connections(phs002455_processed_UMItable, nodes = "Cluster_Refined", group_by = "Skin", print_progress = T)
load(url("https://www.dropbox.com/s/p30mwg5pcj58cc6/network_table.Rdata?dl=1"))
str(network_table)
List of 3
$ Summary :'data.frame': 507 obs. of 6 variables:
..$ Lig_produce : chr [1:507] "CD8" "CD8" "CD8" "CD8" ...
..$ Rec_receive : chr [1:507] "CD8" "CD8" "CD8" "DC" ...
..$ Skin : chr [1:507] "1-Healthy" "2-NonLesional" "3-Disease" "1-Healthy" ...
..$ num_connections : int [1:507] 133 141 141 113 115 115 107 106 110 69 ...
..$ fraction_connections : num [1:507] 0.01154 0.01291 0.01059 0.00981 0.01053 ...
..$ proportion_connections: num [1:507] 0.0968 0.1004 0.0999 0.0822 0.0819 ...
$ full_network :'data.frame': 47540 obs. of 21 variables:
..$ Cluster_Refined : chr [1:47540] "KRT-SP" "KRT-SP" "KRT-SP" "KRT-SP" ...
..$ Ligand : chr [1:47540] "GPI" "HLA-A" "HLA-B" "APP" ...
..$ Cluster_Refined : chr [1:47540] "KRT-SP" "KRT-SP" "KRT-SP" "KRT-SP" ...
..$ Receptor : chr [1:47540] "AMFR" "APLP2" "CANX" "CAV1" ...
..$ Skin : chr [1:47540] "1-Healthy" "1-Healthy" "1-Healthy" "1-Healthy" ...
..$ Ligand_expression : num [1:47540] 142.5 412.4 959.8 58.1 59.9 ...
..$ Receptor_expression : num [1:47540] 36.1 31.3 190.2 81.4 81.4 ...
..$ Connection_product : num [1:47540] 5150 12928 182571 4726 4876 ...
..$ log10_Connection_product : num [1:47540] 3.71 4.11 5.26 3.67 3.69 ...
..$ Connection_ratio : num [1:47540] 3.946 13.157 5.046 0.714 0.736 ...
..$ log10_Connection_ratio : num [1:47540] 0.596 1.119 0.703 -0.146 -0.133 ...
..$ Connection_rank_coarse : num [1:47540] 1913 1361 308 1971 1957 ...
..$ Connection_Z_coarse : num [1:47540] 2.46 2.76 3.64 2.43 2.44 ...
..$ Connection_rank_coarse_grouped: num [1:47540] 684 493 104 712 705 ...
..$ Connection_Z_coarse_grouped : num [1:47540] 2.41 2.7 3.56 2.38 2.39 ...
..$ Connection_rank_fine : num [1:47540] 92 52 4 99 95 44 63 161 5 187 ...
..$ Connection_Z_fine : num [1:47540] 2.99 3.35 4.38 2.96 2.97 ...
..$ Connection_rank_fine_grouped : num [1:47540] 38 23 2 41 40 19 26 62 3 68 ...
..$ Connection_Z_fine_grouped : num [1:47540] 2.85 3.2 4.18 2.82 2.83 ...
..$ Zscores_genes : num [1:47540] 0.0813 -0.3408 -1.5589 0.9897 1.4664 ...
..$ Zscores_genes_grouped : num [1:47540] 0.0538 -0.3484 -1.528 0.9381 1.5023 ...
$ full_network_raw:'data.frame': 351351 obs. of 21 variables:
..$ Cluster_Refined : chr [1:351351] "KRT-SP" "KRT-SP" "KRT-SP" "KRT-SP" ...
..$ Ligand : chr [1:351351] "CALM1" "CALM2" "CALM3" "LIN7C" ...
..$ Cluster_Refined : chr [1:351351] "KRT-SP" "KRT-SP" "KRT-SP" "KRT-SP" ...
..$ Receptor : chr [1:351351] "ABCA1" "ABCA1" "ABCA1" "ABCA1" ...
..$ Skin : chr [1:351351] "1-Healthy" "1-Healthy" "1-Healthy" "1-Healthy" ...
..$ Ligand_expression : num [1:351351] 409.5 608.3 203.5 32.3 54.3 ...
..$ Receptor_expression : num [1:351351] 0 0 0 0 0 0 0 0 0 0 ...
..$ Connection_product : num [1:351351] 0 0 0 0 0 0 0 0 0 0 ...
..$ log10_Connection_product : num [1:351351] 0 0 0 0 0 0 0 0 0 0 ...
..$ Connection_ratio : num [1:351351] 0 0 0 0 0 0 0 0 0 0 ...
..$ log10_Connection_ratio : num [1:351351] 0 0 0 0 0 0 0 0 0 0 ...
..$ Connection_rank_coarse : num [1:351351] 15109 15109 15109 15109 15109 ...
..$ Connection_Z_coarse : num [1:351351] -0.356 -0.356 -0.356 -0.356 -0.356 ...
..$ Connection_rank_coarse_grouped: num [1:351351] 5035 5035 5035 5035 5035 ...
..$ Connection_Z_coarse_grouped : num [1:351351] -0.357 -0.357 -0.357 -0.357 -0.357 ...
..$ Connection_rank_fine : num [1:351351] 1145 1145 1145 1145 1145 ...
..$ Connection_Z_fine : num [1:351351] -0.327 -0.327 -0.327 -0.327 -0.327 ...
..$ Connection_rank_fine_grouped : num [1:351351] 382 382 382 382 382 382 382 382 382 382 ...
..$ Connection_Z_fine_grouped : num [1:351351] -0.329 -0.329 -0.329 -0.329 -0.329 ...
..$ Zscores_genes : num [1:351351] -0.288 -0.288 -0.288 -0.263 -0.191 ...
..$ Zscores_genes_grouped : num [1:351351] -0.288 -0.288 -0.288 -0.263 -0.19 ...
# The lists contains 3 entries.
# Summary provides a high level view of the number of connections between cell types.
# full_network is the same as full_network_raw, except that full_network_raw will preserve entries where either the ligand or receptor is 0 in expression.
head(network_table$full_network)
Cluster_Refined Ligand Cluster_Refined Receptor Skin Ligand_expression
27 KRT-SP GPI KRT-SP AMFR 1-Healthy 142.54554
28 KRT-SP HLA-A KRT-SP APLP2 1-Healthy 412.41139
46 KRT-SP HLA-B KRT-SP CANX 1-Healthy 959.83108
48 KRT-SP APP KRT-SP CAV1 1-Healthy 58.08089
50 KRT-SP CYR61 KRT-SP CAV1 1-Healthy 59.91347
51 KRT-SP GNAI2 KRT-SP CAV1 1-Healthy 195.81901
Receptor_expression Connection_product log10_Connection_product
27 36.12779 5149.856 3.711795
28 31.34619 12927.527 4.111515
46 190.21150 182570.911 5.261432
48 81.37710 4726.455 3.674535
50 81.37710 4875.584 3.688027
51 81.37710 15935.183 4.202357
Connection_ratio log10_Connection_ratio Connection_rank_coarse
27 3.9455924 0.5961122 1913
28 13.1566661 1.1191459 1361
46 5.0461254 0.7029580 308
48 0.7137253 -0.1464689 1971
50 0.7362448 -0.1329777 1957
51 2.4063160 0.3813527 1231
Connection_Z_coarse Connection_rank_coarse_grouped
27 2.459681 684
28 2.762894 493
46 3.635177 104
48 2.431417 712
50 2.441651 705
51 2.831803 433
Connection_Z_coarse_grouped Connection_rank_fine Connection_Z_fine
27 2.405913 92 2.994012
28 2.703450 52 3.351666
46 3.559407 4 4.380566
48 2.378178 99 2.960674
50 2.388220 95 2.972745
51 2.771069 44 3.432948
Connection_rank_fine_grouped Connection_Z_fine_grouped Zscores_genes
27 38 2.853324 0.08134314
28 23 3.196009 -0.34080106
46 2 4.181846 -1.55893737
48 41 2.821381 0.98973688
50 40 2.832947 1.46642871
51 19 3.273889 0.80125671
Zscores_genes_grouped
27 0.05375211
28 -0.34838912
46 -1.52802185
48 0.93809202
50 1.50234668
51 0.77596199
Each row is defined as a cell type expressing a ligand, a cell type expressing the receptor, as well as the condition (Skin, ie Healthy, Non-lesional, Lesional), the expression values of the ligand and receptor, as well as the log10_Connection_product, which is a log10 transformed product of the ligand and expression CPM values.
There are also some basic statistics for these connections (z scores, etc), however these are not used in the manuscript and are beyond the scope of this vignette.
Now we can take this long format table, and convert it into an igraph object. There are several ways this could be done…
One way would be to create 3 separate networks. One for healthy, one for non-lesional, one for lesional. However this has caveats. Because a ligand or receptor may go from OFF to ON, this means the landscape (the nodes and edges of the network) would change between conditions. To avoid this complexity, we opted to create a network that is a superset of all available networks (merge_all = T argument). In this way all possible nodes and edges are represented. Later on, we then query the original network_table in order to determine the edge values associate with each condition.
Please note that these functions write out files by DEFAULT. These files include plots and data. Set your working directory accordingly to track these files location. The prefix argument determines the prefix for these written files. Here we use the log10_Connection_product as the edge values for this network.
network_igraph <- build_rl_network(input = network_table, merge_all = T, value = "log10_Connection_product", prefix = "scitranslmed.abd8995_network")
names(network_igraph)
[1] "igraph_Network" "layout" "clusters"
[4] "clusters_subgraphs" "interactive"
The igraph_network is the full network
The layout is a 2D matrix with node positions
Clusters shows the cluster membership of each node
clusters_subgraphs contains the subgraphs for each cluster
interactive is the web browser based display. Keep in mind this renders in browser and can take a long time to fully render for the full network.
The network contains both a main connected body, as well as some orphan ligands and receptors with no annotated edge to connect them to the main network. For this reason, further clustering is only performed on the main graph (subset = 1, subset_on = “connected” arguments)
network_igraph_C1_analyzed <- analyze_rl_network(network_igraph, subset = 1, subset_on = "connected", prefix = "scitranslmed.abd8995_network_C1", cluster_type = "louvain", merge_singles = F)
names(network_igraph_C1_analyzed)
[1] "Degree" "Node_degree"
[3] "Node_betweeness" "Node_authority"
[5] "Edge_degree" "Edge_betweeness"
[7] "Edge_hub" "Crossing"
[9] "Clusters_Results" "Clusters_individual"
[11] "Clusters_betweeness" "Clusters_edge_hub"
[13] "Clusters_node_authority" "Interactive"
[15] "igraph_Network" "layout"
The output contains statistics related to the clustering results, as well as individual cluster results and the interactive network.
The html files are a great way to view the network. Please note that they are rendered in the browser, so the full network can take a few minutes to display properly. There are drop down menus that allow you to select a particular node of interest, as well as a particular community (cluster) of the network.
input_cell_net <- network_igraph_C1_analyzed$igraph_Network
V(input_cell_net)$membership <- network_igraph_C1_analyzed$Clusters_Results$membership
clusters <- network_igraph_C1_analyzed$Clusters_Results
weights_clusters <- ifelse(crossing(clusters, input_cell_net), 1, 10)
set.seed(10)
cluster_weighted_layout <- layout_with_fr(input_cell_net, weights = weights_clusters)
edge_colors <- ifelse(crossing(clusters, input_cell_net), "red", "black")
members <- V(input_cell_net)$membership
colopal <- RColorBrewer::brewer.pal(n = 8, name = "Dark2")
f <- colorRampPalette(colopal)
colopal <- f(length(unique(members)))
colopal <- colopal[as.factor(members)]
plot(input_cell_net,
layout = cluster_weighted_layout,
vertex.color = colopal,
vertex.size = 1.5,
vertex.label = "",
vertex.label.cex = 1,
vertex.label.color = "black",
vertex.frame.color=colopal,
edge.width = 0.1,
edge.arrow.size = 0.01,
edge.arrow.width = 0.1,
edge.color = "gray50")